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Abstract 

The present work is devoted to study the numerical simulation for unsteady MHD flow and heat transfer of a couple stress 
fluid over a rotating disk. A similarity transformation is employed to reduce the time dependent system of nonlinear partial 
differential equations (PDEs) to ordinary differential equations (ODEs). The Runge-Kutta method and shooting technique are 
employed for finding the numerical solution of the governing system. The influences of governing parameters viz. 
unsteadiness parameter, couple stress and various physical parameters on velocity, temperature and pressure profiles are 
analyzed graphically and discussed in detail. 
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Introduction 

The flow problems of non-Newtonian fluids are attracted the 
interest of many researchers because of its significance in modern 
technology and industries. The natural and industrial applications 
of such fluids are volcanic lava, extrusion of polymer fluids, drilling 
mud, suspension solutions, cosmetic and food products, solidifica- 
tion of liquid crystals, cooling of metallic plates in a bath, exotic 
lubricants, colloidal and many others. In recent times, among 
various non-Newtonian fluid models, the couple stress fluid model 
has got the special status because of the spin field in the fluid which 
sets up an anti-symmetric stresses, known as couple stresses. The 
basic theory and constitutive equations for couple stress fluid 
initially developed by Stokes [1], is one amongst the polar fluid 
theories which consider couple stresses in addition to the classical 
Gauchy stresses. This type of fluid encloses inflexible and random- 
oriented particles in a viscous medium. Examples of couple stress 
fluids are liquid crystals, colloidal fluids, liquids containing long- 
chain molecules such as polymer suspensions, blood, lubrications 
and additive electro-rheological fluids etc. These fluids are 
considered by many investigators because of the simplest model 
of couple stress fluid amongst others. (See recent attempts [2-4]). 

Magnetohydrodynamics (MHD) occurs when a conducting fluid 
flows in the presence of an electromagnetic field, such that the 
induced electromagnetic force affects the motion of the fluid. 
Some significant technological applications of MHD are, in 
experiments of controlled thermonuclear fusion where a strong 
magnetic field is used to confine rings or columns of hot plasma, to 
generate electricity where liquid metals are driven through a 
magnetic field, etc. Moreover, MHD principles are used in plasma 
accelerators for ion thrusters, for spacecraft propulsion and for 
light-ion-beam powered inertial confinement. Some remarkable 
attempts have been made to study the influence of MHD on 
various flow situations. The influences of MHD stagnation point 



flow was analyzed by Hayat et al. [5], MHD flow in a circular 
magnetic field was considered by Kumar et al. [6], viscous 
dissipation and Joule heating on unsteady mixed convection MHD 
flow was discussed by Osalusi et al. [7] . Investigation of fractional 
MHD Oldroyd-B fluid over an oscillating plate was the subject of 
Jamil et al. [8], and the influence of MHD in the porous medium 
and heat transfer analysis of couple stress fluid between two 
parallel plates with variable viscosity was investigated by Sreenadh 
et al. [9]. Numerical methods to obtain the solutions of some 
complicated MHD flow problems were presented in [10-12]. 

Rotating-disk flow of non-Newtonian fluids is intriguingly 
growing interest, as it is momentous in the study of engineering 
flows on rotating machineries, spin-coating, centrifugal pumps, 
rotational viscometers, computer storage devices and some 
aerodynamic related problems in fluid mechanics. Taking into 
account both industrial and technological applications, there have 
been a number of studies presented in literature on steady flow of 
various types of fluid over a rotating disks and sphere with and 
without considering magnetic effect [13-25]. It is worth mention- 
ing that the real world applications encounter the unsteadiness in 
the flow and temperature fields of almost all processes. Osalusi 
et al. [26] presented the eflects of Ohmic heating and viscous 
dissipation on unsteady MHD and slip flow over a porous rotating 
disk. The series solutions for unsteady laminar MHD flow near 
forward stagnation point of an impulsively rotating and translating 
sphere in presence of buoyancy forces are found by Dinarvand 
et al. [27]. Ghamkha et al. [28] studied the unsteady MHD flow of 
mixed convection fluid with heat and mass transfer. Ram et al. 
[29] analyzed the influence of phase difference between highly 
oscillating magnetic field and magnetization on the unsteady ferro 
fluid flow due to a rotating disk. Nadeem et al. [30] worked out on 
the analytical treatment of unsteady mixed convection MHD flow 
on a rotating cone in a rotating frame. 
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Figure 1. Coordinates system and flow configuration. 

doi:10.1371/journal.pone.0095423.g001 



Takhar et al. [31] examined the unsteady MHD flow of an 
ambient fluid with heat transfer over a rotating disk, present study 
aims at extending aforesaid work to incorporate the couple stress 
effects with magnetic field and heat transfer analysis is also taken 
into account. However, to the best of our knowledge, no attempt is 
available in the literature regarding unsteady MHD flow of couple 
stress fluid over a rotating disk. In this note, an investigation has 
been carried out to analyze the unsteady flow and heat transfer of 
couple stress fluid with the effect of transverse magnetic field on an 
infinite rotating disk along z-axis. The governing time dependent 
momentum and energy partial differential equations are trans- 
formed into a system of ordinary differential equations by using 
suitable similarity transformations. These ODEs are then solved 
numerically by employing a Runge-Kutta technique coupled with 
shooting method. The variations in the flow field, heat transfer and 
pressure of the fluid due to pertinent parameters are obtained 
through graphs and discussed thoroughly. 



The governing partial diflerential equations are 
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Problem Formulation 

Consider the three dimensional, unsteady, laminar, incom- 
pressible, MHD flow and heat transfer of couple stress fluid over a 
rotating disk which rotates about the z-axis with angular velocity Q 
in a cylindrical polar coordinate system r, 0 and z. The r-axis is in 
the direction of radius of disk, and the applied uniform magnetic 
field ^0 is in the direction of z-axis, normal to the surface of disk. 
As the flow is axisymmetric therefore radial, tangential and axial 
components of velocity w, v, w, temperature T and pressure P are 
independent of 6. The induced magnetic field can be neglected in 
comparison with the applied magnetic field by considering the 
magnetic Reynolds number RCm much less than 1. Fig. 1 
represents the flow description and geometrical configuration of 
the physical model. 
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Figure 2. The radial velocity field for M = 1.5, (^ = 0.8. 

doi:1 0.1 371 /journal. pone.0095423.g002 



where 



K= - 

1 d 

r dr 



The boundary conditions are 



w(r,z,0 = 0, v(r,z,0 = rQ, w(r,z,0 = 0, 
Uzz(r,z,t) = Vz(r,z,t) = 0, 



7V = 



(5) 



T(r,z,t) = T^, asz = 0, 



w(r,z,0 = 0, w(r,z,0 = 0, ^^(^^,0 = ^z(^^,0 = 0, 
T(r,z,t) = Ta3, as z^oo. 



(6) 



(7) 



Here, M, v and w are the components of velocity along the r, 6 
and z directions respectively, P and T are the pressure and 
temperature, respectively, p, v and a are density, viscosity and 
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Figure 3. The radial velocity field for M = l, a = 0.2, 

doi:1 0.1 371 /journal. pone.0095423.g003 
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Figure 4. The radial velocity field for « = 0.5, C = 0.8. 

doi:1 0.1 371 /journal. pone.0095423.g004 

thermal diffusivity, respectively. Bq is the applied magnetic field, t 
is the fDced time and K is the couple stress parameter respectively. 
Using similarity transformations 
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The equation of continuity yields the relation 
2/ + /z' = 0, 



(8) 



(9) 




Figure 5. The azimuthal velocity field for M = 1.5, C = 0.8. 

doi:1 0.1 371 /journal. pone.0095423.g005 
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Figure 6. The azimuthal velocity field for M = l, a = 0.2, 

doi:1 0.1 371 /journal. pone.0095423.g006 

And Eqs. (2) — (5), after using Eqs. (8) — (9) reduce to 
Ch'"" -h"' + a (hh" - y + 2^^ + Mh'\ -{h'+^h"^=0, (10) 



F + ^h""-h" + ahh'-^h'- ^ =0, (12) 



Cg"" - g" + a{hg' - gA' + Mg) - (g + I g') = 0, (11) 



(13) 



with boundary conditions 
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Figure 7. The azimuthal velocity field for a = 0.5, C = 0.8. 

doi:10.1371/journal.pone.0095423.g007 
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Figure 8. The axial velocity field for M = 1.5, C = 0.8. 

doi:10.1371/journal.pone.0095423.g008 



h = h' = 0, h'" = 0, g = (l) = \, g' = 0, P = 0, at f/ = 0, 
/z = 0, h" = 0, g = 0, g' = 0, at rj^ CO. 



(14) 



approximate similarity. Therefore C = valid approxima- 

tion. And the solutions obtained here are approximations to 
solutions of the original nonlinear PDEs. 

The shear stress of the fluid in the radial and axial directions can 
be found as 



Here h\ h and g are the radial, axial and tangential components 
of the dimensionless velocity, (j) and P are the dimensionless 
temperature and pressure of the flow, a is the dimensionless 
unsteadiness parameter which measures the decrease or increase 
in the angular velocity Q, fi is the coefficient of viscosity, Pr is the 
Prandtl number. Ha is the Hartmann number, M is the magnetic 
parameter, RCr is the Reynolds number with respect to radius 
respectively. The dimensionless couple stress parameter C depends 
on time t, the Eqs.(lO) — (13) are not pure similarity equations, to 
avoid ambiguity, t has been fixed which led the above equation to 
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Figure 9. The axial velocity field for M = l, a = 0.2. 

doi:1 0.1 371 /journal. pone.0095423.g009 
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Figure 10. The axial velocity field for fl = 0.5, C = 0.8. 

doi:1 0.1 371 /journal.pone.0095423.g01 0 



Hence the local skin friction coefficients are given by 



C/Re,2 = a4^-r(0)-^r"(0) ), 
C,Re,2 = a-i(g'(0) + &"'(0)). 



(17) 



The heat transfer coefficient in terms of the Nusselt Nu number 
can be expressed as 



Nu = 



dT 



Tw Trr 



-RerL"2f(0), (18) 



We use the result xqz with a radial strip integration to find the 
total torque required to turn a disk of radius b. 



2 = 271 I zei ydr=^^^{m + ^g"'m, (19) 
2t{tv)2 




Figure 11. The temperature gradient for M = \, Pv = \, il, = \, 

doi:1 0.1 371 /journal. pone.0095423.g011 



PLOS ONE I www.plosone.org 



7 



May 2014 | Volume 9 | Issue 5 | e95423 



Flow of Couple Stress Fluid over a Rotating Disk 




Figure 12. The temperature gradient for M = l, a = 0.1, C = 

doi:1 0.1 371 /journal.pone.0095423.g01 2 



The dimensionless moment coefficient for a disk of radius b can 
be defined as / ? /2 



AO 11 



= - ^-f- = - 2 71 Rer2 (g'(0) + &'"(0)) . (20) 



h"-a(^h"-^ +2g^ + Mh^ + =0, (21) 



The system of Eqs. (10) — (12) have been modeled the unsteady / ^ \ 

flow of a couple stress fluid over a rotating disk. The case when 

g"-a(hg'-gh' + Mg)+[g+ jg'j =0, (22) 

there is no body couples or rotation in the fluid, we set the couple 
stress parameter C equal to zero in the governing system, which 
results in the following form 



2.5 




Figure 13. The pressure gradient for M = 1.5, C = 0.6. 

doi:1 0.1 371 /journal.pone.0095423.g01 3 
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V 



Figure 14. The pressure gradient for M = l, a = 0.2. 

doi:10.1371/journal.pone.0095423.g014 



P'-h" + ahh'-^h'- - =0. (23) 

It can be noticed that the Eqs. (21) — (23) represents the model 
of unsteady flow of an electrically conducting viscous fluid which 
was considered and studied by Takhar et. al. [31]. 



Numerical Procedure 

The resulting highly nonlinear system of ODEs (10) — (13) with 
boundary conditions (14) are solved numerically using fourth- 
order Runge-Kutta integration scheme with shooting technique 
for fixed time /q. In shooting method, the resulting system of 
nonlinear ordinary differential equations is transformed to a set of 
simultaneous first order differential equations. A systematic 
guessing of initial conditions and numerical integration as an 
initial value problem has been done until the desired degree of 
accuracy has obtained. We have relied on the computational 
software MATHEMATICA package to solve the considered 




Figure 15. The pressure gradient for a = 0.5, C = 0.8. 

doi:l 0.1 371 /journal.pone.0095423.g01 5 
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problem. Influences of all relevant parameters are presented 
graphically. 

Results and Discussion 

The numerical solution of the problem has been computed via 
fourth order Runge-Kutta method based on shooting technique. 
The effect of various physical parameters namely the unsteadiness 
parameters a, magnetic parameter M, couple stress parameter C 
and Prandtl number Pr have been examined to measure the 
variations in the angular velocity, in the velocity field, in the 
temperature gradient and in the pressure distribution of the fluid 
(see Figs. 2-15). The changes in radial velocity / for the variations 
in the parameters a, M and C are displayed in Figs. 2-4. It can be 
seen from Fig. 2 that the unsteadiness parameter a decreases the 
radial velocity /, the stress due to the micro rotation particles in 
the fluid also weakens the motion of the couple stress fluid as 
increase in the couple stress parameter C decreases the radial 
velocity (see Fig. 3), whereas Fig. 4 shows that the increase in 
magnetic parameter M results in the decline of velocity function in 
radial direction which indicates that the inhibiting Lorentz forces 
opposes the flow in this direction. 

The effects of tangential velocity g as a function of rj for a, M 
and C are shown in Figs. 5-7. The decline in the tangential velocity 
component can be observed with the increase in unsteadiness 
parameter a as shown in Fig. 5. An increase in the tangential flow 
over the disk is observed with the increment in the couple stress 
parameter C increasing the boundary wall thickness (see Fig. 6). 
Since the magnetic force is in the direction normal to the disk, it 
opposes the tangential flow velocity, therefore the tangential 
velocity is lowered everywhere. Thus the effect of the Lorentz 
forces leads to decelerate the azimuthal velocity with increase in 
M. (see Fig. 7). 

Figs. 8-10 shows the effect of aforesaid parameters on axial 
velocity A as a function of f]. The axial velocity field behaves in a 
similar way to the radial velocity field. Flow field is reduced in the 
axial direction with the enhancement in the unsteadiness 
parameter a (as shown in Fig. 8), Fig. 9 shows the effects of 
increasing couple stress parameter which brings unfavorable 
change in the axial velocity component. The effect of magnetic 
parameter can be observed in Fig. 1 0 which displays the fall in the 
axial velocity component of couple stress fluid flow with the 
increase in magnetic parameter M. 

In Figs. 11-12, the effect of temperature gradient (/> as a function 
of f] on the different values of Pr and a are observed. Increasing a 
results in decline in the temperature gradient (Fig. 11), while the 
similar behavior is observed with increasing Pr (Fig. 12). 

Figs. 13-15 describes the outcome of the pressure distribution 
against a, M and C- These figures depicts the inverse relation 
between pressure distribution and the unsteadiness parameter, 
magnetic parameter and couple stress parameter, indicating an 
increase in above mentioned parameters led the decreasing 
behavior of the pressure distribution. 

Limitations of the Study and Future 
Recommendations 

It is imperative to enlighten diverse limitations of this research. 
This discussion aims at assisting readers to comprehend this study, 
which ultimately, provide an opportunity to extend the current 
research. The following assumptions and limitations are consid- 
ered. 

• Flow is unsteady, incompressible and laminar. 



• Flow is three dimensional and rotating over a disk about the 
z— axis. 

• A uniform magnetic field is applied in the direction normal to 
the surface of the disk which is assumed unchanging with a 
small magnetic Reynolds number RCm « 1 . 

• It is assumed that the effect of dissipation in the energy 
equation is negligible. 

In this study we attempt to solve governing equations 
numerically for a fixed value of (i-e., at a fixed time ) since 

C(^ — • interesting to note that for large time the higher 

derivatives in equations (10) — (11) are vanished as C tends to zero. 
This leads the case of degenerate limit and difficult to solve with 
the presented technique, so a different approach should be sought 
such as singular perturbation theory to solve the problem for the 
large-time limit. It is worth observing that, on considering the 
aforesaid limitations for future research, it will turn out to be a 
demanding task for the researchers regardless of their particular 
interests in this area. For the model like couple stress fluid, an 
extensive range of problems can be defined for both theoretical 
and experimental aspects, i.e. flow over an off centered disk, 
between two coaxial disk, due to rotating sphere, over an inclined 
disk etc. So, fascinated by the diverse applications of non- 
Newtonian fluids, the present work can be extended to other non- 
Newtonian fluids model by employing a variety of numerical 
methods. 

Conclusions 

In this paper, the response of unsteady MHD flow of couple 
stress fluid over a rotating disk is discussed. The nonlinear 
differential equations have been solved numerically using Runge- 
Kutta method and shooting technique. The numerical formulas 
for radial and tangential shear stress, heat transfer coefficient in 
term of Nusselt number and dimensionless moment coefficient are 
presented. The influence of various physical parameters and 
couple stress parameter on the velocity field, temperature and 
pressure distributions are all examined numerically and discussed 
graphically. The following results can be drawn: 

• The unsteadiness parameter decelerates the flow and temper- 
ature field over the rotating disk. 

• The presence of couple stresses in the fluid enhances the fluid 
motion in tangential direction and weakens in radial and axial 
directions. 

• The magnetic field opposes the fluid velocity in all directions 
i.e. radial, tangential and axial. 

• There is a significant reduction in the temperature field due to 
an increase in Prandtl number and unsteadiness parameter. 

• Pressure gradient decline while unsteadiness parameter, 
magnetic parameter and couple stress parameter increases. 
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